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Abstract 

The relativistic quasiparticle random phase approximation (RQRPA) is formulated in the canon- 
ical single-nucleon basis of the relativistic Hartree-Bogoliubov (RHB) model. For the interaction 
in the particle-hole channel effective Lagrangians with nonlinear meson self-interactions are used, 
and pairing correlations are described by the pairing part of the finite range Gogny interaction. 
The RQRPA configuration space includes the Dirac sea of negative energy states. Both in the 
particle-hole and particle-particle channels, the same interactions are used in the RHB calculation 
of the ground state and in the matrix equations of the RQRPA. The RHB+RQRPA approach is 
tested in the example of multipole excitations of neutron rich oxygen isotopes. The RQRPA is 
applied in the analysis of the evolution of the low-lying isovector dipole strength in Sn isotopes 
and N=82 isotones. 

PACS numbers: 21.30.Fe, 21.60.Jz, 24.30.Cz, 24.30.Gd 
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I. INTRODUCTION 



The multipole response of unstable nuclei far from the line of /3-stability presents a very 
active field of research, both experimental and theoretical. These nuclei are characterized 
by unique structure properties: the weak binding of the outermost nucleons and the effects 
of the coupling between bound states and the particle continuum. On the neutron rich side, 
in particular, the modification of the effective nuclear potential leads to the formation of 
nuclei with very diffuse neutron densities, to the occurrence of the neutron skin and halo 
structures. These phenomena will also affect the multipole response of unstable nuclei, 
in particular the electric dipole and quadrupole excitations, and new modes of excitations 
might arise in nuclei near the drip line. 

A quantitative description of ground-states and properties of excited states in nuclei 
characterized by the closeness of the Fermi surface to the particle continuum, necessitates a 
unified description of mean-field and pairing correlations, as for example in the framework 
of the Hartree-Fock-Bogoliubov (HFB) theory. In order to describe transitions to low-lying 
excited states in weakly bound nuclei, in particular, the two-quasiparticle configuration space 
must include states with both nucleons in the discrete bound levels, states with one nucleon 
in a bound level and one nucleon in the continuum, and also states with both nucleons in 
the continuum. This cannot be accomplished in the framework of the BCS approximation, 
since the BCS scheme does not provide a correct description of the scattering of nucleonic 
pairs from bound states to the positive energy particle continuum. Collective low-lying 
excited states in weakly bound nuclei are best described by the quasiparticle random phase 
approximation (QRPA) based on the HFB framework. The HFB based QRPA has been 
investigated in a number of recent theoretical studies. In Ref. |fl] a fully self-consistent 
QRPA has been formulated in the HFB canonical single-particle basis. The Hartree-Fock- 
Bogoliubov formalism in coordinate state representation has also been used as a basis for 
the continuum linear response theory 

3,3- In Ref. Q the HFB energy functional has been 
used to derive the continuum QRPA response function in coordinate space. HFB based 
continuum QRPA calculations have been performed for the low-lying excited states and 
giant resonances, as well as for the j3 decay rates in neutron rich nuclei. 

In this work we formulate the relativistic QRPA in the canonical single-nucleon basis of 
the relativistic Hartree-Bogoliubov (RHB) model. The RHB model is based on the rela- 



2 



tivistic mean-field theory and on the Hartree-Fock-Bogoliubov framework. It has been very 
successfully applied in the description of a variety of nuclear structure phenomena, not only 
in nuclei along the valley of /^-stability, but also in exotic nuclei with extreme isospin values 
and close to the particle drip lines. Another relativistic model, the relativistic random phase 
approximation (RRPA), has been recently employed in quantitative analyses of collective 
excitations in nuclei. Two points are essential for the successful application of the RRPA in 
the description of dynamical properties of finite nuclei: (i) the use of effective Lagrangians 
with nonlinear self-interaction terms, and (ii) the fully consistent treatment of the Dirac sea 
of negative energy states. 

The RRPA with nonlinear meson interaction terms, and with a configuration space that 
includes the Dirac sea of negative-energy state, has been very successfully employed in 
studies of nuclear compressional modes J5, , of multipole giant resonances and of low- 
lying collective states in spherical nuclei |S | , of the evolution of the low- lying isovector dipole 
response in nuclei with a large neutron excess and of toroidal dipole resonances U|. 

In Sec. iniwe present the formalism and formulate the matrix equations of the relativistic 
QRPA in the canonical basis of the relativistic Hartree-Bogoliubov (RHB) framework for 
spherical even-even nuclei. In Sec. IHII the RHB+RQRPA approach is tested in the example 
of the isoscalar monopole, isovector dipole and isoscalar quadrupole response of 22 0, and 
the results are compared with recent non-relativistic QRPA calculations of the multipole 
response of neutron rich oxygen isotopes. In Sec. IIVI the RQRPA framework is applied in 
the analysis of the evolution of the low-lying isovector dipole strength in Sn isotopes and 
N=82 isotones. The results are compared with recent experimental data. Section M contains 
the summary and the conclusions. 



II. THE RELATIVISTIC QUASIPARTICLE RANDOM PHASE APPROXIMA- 
TION 

In this section the matrix equations of the relativistic quasiparticle random phase approxi- 
mation (RQRPA) are formulated in the canonical basis of the relativistic Hartree-Bogoliubov 
(RHB) framework for spherical even-even nuclei. 
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A. The relativistic mean-field Lagrangian and the equations of motion 



The nuclear matter equation of state and detailed properties of finite nuclei have been very 
successfully described by relativistic mean- field (RMF) models [LI HQ- In this framework 



the nucleus is described as a system of Dirac nucleons that interact in a relativistic covariant 
manner by meson exchange. In particular, the isoscalar scalar cr-meson, the isoscalar vector 
cj-meson, and the isovector vector p-meson build the minimal set of meson fields that is 
necessary for a quantitative description of bulk and single-particle nuclear properties. The 
model is defined by the Lagrangian density 

£ = Cn + £m + Cint- (1) 

Cn denotes the Lagrangian of the free nucleon 

C N = tfj(tYd fl -m)^, (2) 

where m is the bare nucleon mass and ip denotes the Dirac spinor. C m is the Lagrangian of 
the free meson fields and the electromagnetic field 

-\R, V R^ + - \F, U F^, (3) 

with the corresponding masses m a , m p , and fl pu , R pu , F pu are field tensors 

R pv = - dvPfj, (4) 
Fy,v = d^Aj, — d v A p . 

The model Lagrangian density contains also the interaction terms 

C mt = -#>V> - - - jr^A^. (5) 

The vertices read 

T a = g a , r^ = ^ 7 M , ?$ = 9 P ?r, r™ = e±^-j», (6) 

with the coupling constants g a , g^, g P and e. This simple linear model, however, does 
not provide a quantitative description of complex nuclear systems. An effective density 
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dependence has been introduced by replacing the quadratic cr-potential -m 2 a 2 with a 

quartic potential U(a) = -m 2 a 2 + t^ " 3 + ~j~°" 4 - This potential includes the nonlinear a 
self- interactions with two additional parameters g<i and g%. It has been shown that one can 
describe the properties of nuclear matter and finite nuclei with high accuracy using density 
dependent coupling constants g m (p), instead of nonlinear a self- interaction 16]. 

From the model Lagrangian density the classical variation principle leads to the equations 
of motion. The time- dependent Dirac equation for the nucleon reads 

[Y(td^ + V^) + m + S]^ = 0. (7) 

If one neglects retardation effects for the meson fields, a self-consistent solution is obtained 
when the time- dependent mean-field potentials 

S(r,t) = g a cr(r,t) , 

V^{r,t) = gjv^r, t) + gpTp^r, t) + eA^r, t) ~ ^ , (8) 
are calculated at each step in time from the solution of the stationary Klein-Gordon equations 

- A0 m + u'{<f> m ) = ± Ur m ^) , (9) 



where the upper sign holds for vector fields and the lower sign for the scalar field. The index 
m denotes mesons and the photon, i.e. m = a, u^, p 5 ", A 11 . This approximation is justified 
by the large meson masses. The corresponding meson exchange forces are of short range 
and therefore retardation effects can be neglected. 

In practical applications to nuclear matter and finite nuclei, the relativistic models are 
used in the no-sea approximation: the Dirac sea of states with negative energies does not 
contribute to the densities and currents. For a nucleus with A nucleons 

A 

(i>T m ^) = &(r, t)T m ^(r, t) , (10) 
i=i 

where the summation is performed only over the occupied orbits in the Fermi sea of positive 
energy states. The set of coupled equations (JJJ) and Q define the relativistic mean field 
(RMF) model. In the stationary case they reduce to a nonlinear eigenvalue problem, and 
in the time-dependent case they describe the nonlinear propagation of the Dirac spinors in 



The mean-field approximation represents the lowest order of the quantum field theory: 
the meson field operators are replaced by their expectation values. The A nucleons, described 
by a Slater determinant |<E>) of single-particle spinors, move independently in the classical 
meson fields. The couplings of the meson fields to the nucleon are adjusted to reproduce the 
properties of nuclear matter and finite nuclei. The a-meson approximates a large attractive 
scalar field, that is produced by very complicated microscopic processes, such as uncorrelated 
and correlated pion-exchange. The cu-meson describes the short range repulsion between the 
nucleons, and the p-meson carries the isospin quantum number. The latter is required by 
the large empirical asymmetry potential in finite nuclear systems. The basic ingredient 
of the microscopic nuclear force is the pion. In relativistic mean-field models it does not 
contribute on the Hartree level because of parity conservation. The pion field has been 
included in the relativistic Hartree-Fock model. However, the resulting equations of motion 
are rather complicated and this model has been rarely used. Many effects that go beyond 
the mean-field level are apparently neglected in the RMF model. Among them are the Fock 
terms, the vacuum polarization effects and the short range Brueckner-type correlations. The 
experimental data to which the meson-nucleon couplings are adjusted, however, contain all 
these effects and much more. It follows that these effects are not completely neglected. On 
the contrary, they are taken into account in an effective way. The concept behind the RMF 
model is therefore equivalent to that of the density functional theory, which is widely used in 
solid state physics, molecular physics, chemistry and also in non-relativistic nuclear physics. 
The RMF model represent the covariant form of this method. 



B. Covariant density functional theory 

The equations of motion of the relativistic mean field model can also be derived starting 
from a density functional. From the energy-momentum tensor one writes the total energy 
of the nuclear system 



+ 



i=i 



/ 



: {Vaf + U{a) 
(Vp) 2 + mV 



d 3 r - 



(Vw) 2 + mlu 2 ' 



d 3 r 



1 



d 3 r~ - I (VA) 2 d 3 r 
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+ 



9<rPsCr + gj^ + gpjpp* + ej Cfl A^ 



d 3 r. 



(11) 



By using the definition of the relativistic single nucleon density matrix 



/)(r,r',t)=^|^(r,t))(^(r',t)| 



(12) 



the total energy can be written as a functional of the density matrix p and the meson fields 



E 



RMF 



[/V 



Tr [(ocp + (3m) p] 



± 



2(V0 m ) 2 H 
-Tr [{Y m (j) m )p} . 



U{<(> r , 



d 3 r 



(13) 



The trace operation involves a sum over the Dirac indices and an integral in coordinate 
space. The index m is used as generic notation for all mesons and the photon. From the 
classical time-dependent variational principle 

rt 2 



5 [ 2 dt{(1>\id t \<l>)-E{p,<j )m ]} = 

Jt! 



(14) 



the equations of motion (JJJ) and © are obtained. The equation of motion for the density 
matrix reads 



id t p = h(p,(f) m ),p 



(15) 



The single particle Hamiltonian h is the functional derivative of the energy with respect to 
the single particle density matrix p 



5E 



(16) 



C. Pairing correlations and the Relativistic Hartree-Bogoliubov theory 



The inclusion of pairing correlations is essential for a quantitative description open shell 
nuclei. In Ref. [l8( a fully microscopic derivation of the relativistic Hartree-Bogoliubov 
theory has been developed. Using the Gorkov factorization technique, it has been shown 
that the pairing interaction results from the one-meson exchange (a-, u- and p-mesons). 
In practice, however, it turns out that the pairing correlations calculated in this way, with 
coupling constants taken from the standard parameter sets of the RMF model, are too 
strong. The repulsion produced by the exchange of vector mesons at short distances results 



in a pairing gap at the Fermi surface that is by a factor three too large. However, as has been 
argued in many applications of the Hartree-Fock-Bogoliubov theory, there is no real reason 
to use the same effective forces in both the particle-hole and particle-particle channels. 

Pairing correlations can be easily included in the framework of density functional theory, 
by using a generalized Slater determinant |$) of the Hartree-Bogoliubov type. The ground 
state of a nucleus |$) is represented as the vacuum with respect to independent quasiparticle 
operators 

at = Y,U lk ct + V lk c l} (17) 
i 

where U^, are the Hartree-Bogoliubov coefficients. They determine the Hermitian single 
particle density matrix 

p = V*V T , (18) 

and the antisymmetric pairing tensor 

k = V*U t . (19) 

The energy functional depends not only on the density matrix p and the meson fields m , 
but in addition also on the pairing tensor. It has the form 

E[p, k, 4> m \ = E RMF [p, 4> m } + E pair [k], (20) 

where E RMF [p, is the RMF- functional defined in Eq. (fT3*j) . The pairing energy E pair [k] 
is given by 

E pair [k} = ^Ti[k*V pp k}. (21) 

V pp is a general two-body pairing interaction. Finally, the total energy can be written as a 
functional of the generalized density matrix 




K = | ' | , (22) 

which obeys the equation of motion 

idtK= [H(H),K]. (23) 

The generalized Hamiltonian 7i is a functional derivative of the energy with respect to the 
generalized density 

/ 



-H 5E 
6TZ 



hp — m — \ A 

-A* -h D + m + \ 
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(24) 



It contains two average potentials: the self-consistent mean field ho, which encloses all 
the long range particle- hole (ph) correlations, and the pairing field A, which includes the 
particle-particle (pp) correlations. The single particle potential ho results from the variation 
of the energy functional with respect to the Hermitian density matrix p 

h D = f p , (25) 

and the pairing field is obtained from the variation of the energy functional with respect to 
the pairing tensor 

* ■ f ■ 

The pairing field is an integral operator with the kernel 

Aa.(r, r') = \ £ C(r> r')^(r, r'), (27) 

c,d 

where a, b, c, d denote quantum numbers that specify the Dirac indices of the spinors, and 
Vabcd( v -i r are the matrix elements of a general two-body pairing interaction. 

The stationary limit of Eq. (J2*3|) describes the ground state of an open-shell nucleus I20I 



V 





(28) 



2l|. It is determined by the solutions of the Hartree-Bogoliubov equations 

' f ho — m — A A 

—A* -h D + m + \ 

The chemical potential A is determined by the particle number subsidiary condition in order 
that the expectation value of the particle number operator in the ground state equals the 
number of nucleons. The column vectors denote the quasiparticle wave functions, and Ek 
are the quasiparticle energies. The dimension of the RHB matrix equation is two times the 
dimension of the corresponding Dirac equation. For each eigenvector ([/&, Vk) with positive 
quasiparticle energy Ek > 0, there exists an eigenvector (V^,[/£) with quasiparticle energy 
— E^. Since the baryon quasiparticle operators satisfy fermion commutation relations, the 
levels Ek and — Ek cannot be occupied simultaneously. For the solution that corresponds to 
a ground state of a nucleus with even particle number, one usually chooses the eigenvectors 
with positive eigenvalues Ek. 

The RHB equations are solved self-consistently, with potentials determined in the mean- 
field approximation from solutions of static Klein-Gordon equations 



-A + ml 



a(r) = -g a p 8 {v) - g 2 a 2 (r) - g 3 a 3 (r) (29) 
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A + mi u°(t) = g uPv {r) (30) 



-A + m 2 p 



p u (r) = g pP3 (r) (31) 
—AA°(r) = ep p (r) (32) 

for the (T-meson, the tu-meson, the p-meson and the photon field, respectively. Because 
of charge conservation, only the 3-rd component of the isovector p-meson contributes. In 
the ground-state solution for an even-even nucleus there are no currents (time reversal 
invariance) and the spatial components u>, p 3 , A of the vector fields vanish. In nuclei 
with an odd number of protons or neutrons time reversal symmetry is broken, and the 
resulting spatial components of the meson fields play an essential role in the description 
of magnetic moments and of moments of inertia in rotating nuclei. The equation for the 
isoscalar scalar u-meson field contains nonlinear terms. The inclusion of nonlinear meson self- 
interaction terms in meson-exchange RMF models is absolutely necessary for a quantitative 
description of ground-state properties of spherical and deformed nuclei Q|. The source 
terms in equations (J29j) to (J32j) are sums of bilinear products of baryon amplitudes 

Ps(r) = E^M/VUr) (33) 

fc>0 

Mr) = EtoW (34) 

fc>0 

P 3 (r) = £^( r h^(r) (35) 

fc>0 

P em (r) = £^(r)^V*(r), (36) 

where Y^k>o * s a shorthand notation for the no-sea approximation. The self-consistent 
solution of the Dirac-Hartree-Bogoliubov integro-differential equations and Klein-Gordon 
equations for the meson fields determines the ground state of a nucleus. In the present im- 
plementation of the RHB model the coupled system of equations is solved by expanding the 
nucleon spinors £4(r) and Vfc(r), and the meson fields in the spherical harmonic oscillator 



basis 
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D. The Relativistic Quasi- Particle Random Phase Approximation 

In this section we will derive the RQRPA from the time-dependent RHB model in the 
limit of small amplitude oscillations. The generalized density matrix TZ and the fields m = 

10 



er, u^, p^, A^- have been considered as independent variables related only by the equations of 
motion. One can use the Klein-Gordon equations to eliminate the meson degrees of freedom, 
but this is only possible in the small amplitude limit. The time-dependent meson field can 
be written as 

m = tf+# m , (37) 

where (f)}$ is the meson field that corresponds to the stationary ground state, and 5<p m 
is a small variation of the meson field around the stationary state solution. In the linear 
approximation the corresponding Klein-Gordon equation reads 

-A + U"(^)] 50 m (r) = ±g m 5p m (r), (38) 

where Sp m (r) are the various densities and currents (see Eq. (fTU|) ). If there are no nonlinear 
meson self-interaction terms, U"{<p^) = m^. The propagator G m (r, r') can be obtained 
analytically and it has the Yukawa form. In the case of nonlinear meson self-interaction 
terms U"((f)$) depends on the field (j>$, and an analytical solution is no longer possible. 
The propagator G m (r,r') has to be calculated numerically (for details see Ref. [23| ) . In 
both cases we find a linear relation between 5<f) m and Sp m 

50 m (r) = ±g m J rf 3 r'G m (r, r')8p m (r'). (39) 

The generalized Hamiltonian Ti, can now be expressed as a functional of the generalized 
density 1Z only. In the linear approximation the generalized density matrix is expanded 

n = TZo + 5TZ(t), (40) 

where 1Z is the stationary ground-state generalized density. Since 7l(t) is a projector at all 
times, in linear order 

TZ 5n + 5KK Q = 511. (41) 
In the quasiparticle basis the matrices TZq and TLq = TKJZq) are diagonal 

n =\ and Ho =\ . (42) 

\ 1 J \ -E n 

From Eq. (J41j) it follows that the matrix 5TZ has the form 

/ 



5K 



\ 




(43) 
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The linearized equation of motion reduces to 



SH 
8TZ 



571,710 



id t TZ = [H , 8TZ] + 
Assuming an oscillatory solution 

6K(t) = J2$K M e i ^ t + h.c 



(44) 



(45) 



the RQRPA equation is obtained 



V 




yv 



yv 



(46) 



For k < k', I < I' the RQRPA matrix elements read 

8 2 E 



Akk'w — (E k + E k ')5 k iS, 



k'V 



and 



Bkk'iv 



8 2 E 



(47) 



8R* kk ,8R u , 8R* kk ,8R* lv 
If the two-body Hamiltonian is density independent the matrices A and B have the simple 



forms 



24] 



B kk >,u> = — ($| [a k ia k , [H, aci>ai]\ |$) . 
Using the representation of the Hamiltonian in the quasiparticle basis 
H = E + J2 H kl' a t a k' + 7 E H kk>iiiai®ti®ii®i 

kk' 4 kk'W 

+ £ ( H kk'ii>®k <4«M + + /i.e.) + ( H kkiw®k®k'®v®i + h.c) 



(48) 



(49) 



kk'W 



we find 



A kk i W - H kl 8 k qi - H kll 8 kU - H kl ,8 k >i + H k , v 8 k i + H kk , u , 

, kk'll' 



Bkh/w — AH^,^. 



(50) 



In the quasiparticle representation the matrix H n is diagonal, i.e. H k j = E k 8 k i. The 
matrices if 22 and if 40 are rather complicated expressions containing the two-body ph- and 
pp-matrix elements and the coefficients U and V (for details see Ref. 241 ]) 
In the more general case of a density dependent Hamiltonian the same expressions can 
be used, but one has to take into account the rearrangement terms originating from the 
variation of the interaction with respect to the density p. 
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E. The Relativistic QRPA in the canonical basis 



The full RQRPA equations are rather complicated, because they require the evaluation of 
the matrix elements H^., lv and H^., u , in the basis of the Hartree-Bogoliubov spinors U k (r) 
and I4(r). It is considerably simpler to solve these equations in the canonical basis, in which 
the relativistic Hartree-Bogoliubov wave functions can be expressed in the form of BCS-like 
wave functions. In this case one needs only the matrix elements V£y K , x of the residual ph- 
interaction, and the matrix elements V^, xx , of the pairing pp-interaction, as well as certain 
combinations of the occupation factors u K , v K . The numerical details are described in the 
Appendix. In the following we use the indices k, A, k' and A' to denote states in the canonical 
basis. We emphasize that the solution of the relativistic quasi-particle RPA equations in 
the canonical basis does not represent an approximation. We obtain a full solution and the 
results do not depend on this special choice of the basis. 

Taking into account the rotational invariance of the nuclear system, the quasiparticle 
pairs can be coupled to good angular momentum and the matrix equations of the relativistic 
quasiparticle random phase approximation (RQRPA) read 



A J B J 

B" J A' J 



J£V,JM 

yVyJM 




J£V,JM 

yVyJM 



(51) 



For each RQRPA energy uj U) X v and Y v denote the corresponding forward- and backward- 
going two-quasiparticle amplitudes, respectively. The coupled RQRPA matrices in the 
canonical basis read 

A J _ rjlHJ), fT U ( J h H ll ^S 4- H n ( J h 

kk'XX' - k\ °k'X' - -H K 'A °kX> - -H K A' °k'X + ^ K '\> OkX 

+U'xx>V K % x (52) 

B'kk'XX' = 2 (£kk'£aA' — ^kk'^XX'Wkk'XX' 

H KK >xA-iy^' +J V^x ■ (53) 

H 11 denotes the one-quasiparticle terms 

Hll = (u K u x - v K v x )h KX - (u K v x + v K u x )A K x , (54) 

i.e. the canonical RHB basis does not diagonalize the Dirac single-nucleon mean-field Hamil- 
tonian ho and the pairing field A. The occupation amplitudes v k of the canonical states are 
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eigenvalues of the density matrix. V ph and V pp are the particle-hole and particle-particle 
residual interactions, respectively. Their matrix elements are multiplied by the pairing fac- 
tors ^ and (, defined below by the occupation amplitudes of the canonical states. The 
relativistic particle-hole interaction V ph is defined by the same effective Lagrangian density 
as the mean- field Dirac single-nucleon Hamiltonian hp. V ph includes the exchange of the 
isoscalar scalar a-meson, the isoscalar vector cu-meson, the isovector vector p-meson, and 
the electromagnetic interaction. The two-body matrix elements include contributions from 
the spatial components of the vector fields. 



A A' 



vtn'Vxx f° r a 5 an d the time-components u°, p°, A if J is even 
for the space-components u>, p, A if J is odd 

Vkk'Vxx' f° r a i an d the time-components u°, p°, A if J is odd 
for the space components u>, p, A if J is even 



with the //-coefficients defined by 



and 



l]^ K , = U K V K > ± V K U K > 



tfr K , = u K u K , =f v K v K > 



The RQRPA configuration space includes the Dirac sea of negative energy states. In 
addition to the configurations built from two-quasiparticle states of positive energy, the 
RQRPA configuration space must also contain pair-configurations formed from the fully 
or partially occupied states of positive energy and the empty negative-energy states from 
the Dirac sea. The inclusion of configurations built from occupied positive-energy states 
and empty negative-energy states is essential for current conservation and the decoupling of 
spurious states 2t\ ■ In recent applications of the relativistic RPA it has been shown that the 
fully consistent inclusion of the Dirac sea of negative energy states in the RRPA configuration 



space is essential 
giant resonances 



or a quantitative comparison with the experimental excitation energies of 



2,m 

It should be emphasized that the present RQRPA model is fully consistent: the same 
interactions, both in the particle-hole and particle-particle channels, are used in the RHB 
equation that determines the canonical quasiparticle basis, and in the RQRPA equation 
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(|51|). In both channels the same strength parameters of the interactions are used in the 
RHB and RQRPA calculations. No additional adjustment of the parameters is needed in 
RQRPA calculations. This is an essential feature of our calculations and it ensures that 
RQRPA amplitudes do not contain spurious components associated with the mixing of the 
nucleon number in the RHB ground state (for + excitations), or with the center-of-mass 
translational motion (for 1~ excitations). 

In the next section we present results of illustrative RQRPA calculations of the multipole 
response in spherical nuclei. For the multipole operator Qx^ the response function R(E) is 
defined 



R(E,J)=Y, B(J, 



, 1 



r/2 



vr (e - u v y + (r/2^ 

where T is the width of the Lorentzian distribution, and 



(55) 



B(J,u„ 



Y,\K'f(K\\Qj\W) 



+ (-iy^'+ J Y^ ( K '\\Qj\\ K )\(u K v K/ + (-l) J w) 



(56) 



In all the examples considered in Sec. 11111 the discrete strength distribution are folded by a 
Lorentzian of width r=lMeV . For the state \J,v), the RQRPA transition density reads 



kk' ^ J 



(57) 



where k and k denote the quantum numbers of the large and small components of the 
Dirac spinors, respectively. f K (r) and g K (r) are the corresponding large and small radial 
components. 



III. ILLUSTRATIVE CALCULATIONS AND TESTS OF THE RQRPA 

Nuclear properties calculated with the RHB+RQRPA model will, of course, crucially 
depend on the choice of the effective RMF Lagrangian in the p/i-channel, as well as on the 
treatment of pairing correlations. The most successful RMF effective interactions are purely 
phenomenological, with parameters adjusted to reproduce the nuclear matter equation of 
state and a set of global properties of spherical closed-shell nuclei. In most applications of 
the RHB model, in particular, we have used the NL3 effective interaction 2^] for the RMF 
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effective Lagrangian. Properties calculated with NL3 indicate that this is probably the best 
nonlinear effective interaction so far, both for nuclei at and away from the line of /3-stability. 
In the pp-channel of the RHB model we have used a phenomenological pairing interaction, 
the pairing part of the Gogny force, 



V pp (l,2) = 2 e -((ri-r 2 )M) ( W{ + BiP a - HiP T - MiP a P T ), (58) 

i=l,2 

with the set D1S j^] for the parameters Wi, ifj and M, (i = 1,2). This force has been 
very carefully adjusted to the pairing properties of finite nuclei all over the periodic table. 
In particular, the basic advantage of the Gogny force is the finite range, which automatically 
guarantees a proper cut-off in momentum space. All RHB+RQRPA calculations presented 
in this work have been performed with the NL3+D1S combination of effective interactions. 

In order to illustrate the RHB+RQRPA approach and to test the numerical implementa- 
tion of the RQRPA equations, in this section we calculate the isoscalar monopole, isovector 
dipole and isoscalar quadrupole response of 22 0. Similar calculations for the neutron-rich 
oxygen isotopes were recently performed by Matsuo in the framework of the non- 

relativistic continuum linear response theory based on the Hartree-Fock-Bogoliubov formal- 
ism in coordinate state representation. The two theoretical frameworks differ, of course, 
both in the physical contents, as well as in the numerical implementation. The results can, 
nevertheless, be compared at least at the qualitative level. In the HFB+QRPA model of 
Refs. 

as 

a Woods-Saxon parameterization is adopted for the single-particle potential, and 
a Skyrme-type density dependent delta force is used for the residual interaction in the ph- 
channel of the QRPA. Since the calculation of the single-particle potential and p/i-interaction 
is not self-consistent, the interaction strength of the residual interaction is renormalized for 
each nucleus in such a way that the dipole response has a zero-energy mode corresponding to 
the spurious center of mass motion. For the pairing interaction, a density-dependent delta 
force is used both in the calculation of the HFB pairing field for the ground state, and in the 
linear response equation for the excitations. The calculation is consistent in the pp-channel. 
The present RHB+RQRPA calculations are fully self-consistent: the same combination of 
effective interactions, NL3 in the p/i-channel and Gogny D1S in the pp-channel, are used 
both in the RHB calculation of the ground state and as RQRPA residual interactions. The 
parameters of the RQRPA residual interactions have exactly the same values as those used 
in the RHB calculation. 
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In the analysis of Refs. |2j, |3( , Matsuo has illustrated the importance of a consistent treat- 
ment of pairing correlations in the HFB+QRPA framework. The residual pairing interaction 
in the QRPA generates pronounced dynamical correlation effects on the responses through 
pair density fluctuations. Moreover, the energy weighted sum rules are only satisfied if the 
pairing interaction is consistently included both in the static HFB and in the dynamical 
linear response. We have verified that the results obtained in the HFB+QRPA framework 
are also reproduced in the RHB+RQRPA calculations. 

In the left panel of Fig. ^ we display the monopole strength function of the neutron 
number operator in 22 O. There should be no response to the number operator since it is a 
conserved quantity, i.e. the Nambu-Goldstone mode associated with the nucleon number 
conservation should have zero excitation energy. The dashed curve (no dynamical pairing) 
represents the strength function obtained when the pairing interaction is included only in 
the RHB calculation of the ground state, but not in the residual interaction of the RQRPA. 
The solid line (zero response) corresponds to the full RHB+RQRPA calculation, with the 
pairing interaction included both in the RHB ground state, and in the RQRPA residual 
interaction. The same result was also obtained in the HFB+QRPA calculation for 24 in 
Ref. Q: the spurious strength of the number operator appears when the pairing interaction 
is included only in the stationary solution for the ground state, i.e. when the dynamical 
QRPA pairing correlations are neglected. 

The isoscalar strength functions of the monopole operator Y^t=x r i in 22 0, shown in the 
right panel of Fig. ^ correspond to three different calculations: a) the RMF+RRPA calcu- 
lation without pairing, b) pairing correlations are included in the RHB calculation of the 
ground state, but not in the RQRPA residual interaction (no dynamical pairing), and c) the 
fully self-consistent RHB+RQRPA calculation. Just as in the case of the number operator, 
by including pairing correlations only in the RHB ground state a strong spurious response is 
generated below 10 MeV. The Nambu-Goldstone mode is found at zero excitation energy (in 
this particular calculation it was located below 0.2 MeV) only when pairing correlations are 
consistently included also in the residual RQRPA interaction. When the result of the full 
RHB+RQRPA is compared with the response calculated without pairing, one notices that, 
as expected, pairing correlations have relatively little influence on the response in the region 
of giant resonances above 20 MeV. A more pronounced effect is found at lower energies. The 
fragmentation of the single peak at ~ 12.5 MeV reflects the broadening of the Fermi surface 
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by the pairing correlations. 

The isovector strength function (J n = 1~) of the dipole operator 

AT Z Z N 

Q ^ = WTz^ rpYlm ~ WTz^ rnYlm m 

p=l n=l 

for 22 is displayed in the left panel of Fig. [21 In this example we also compare the results of 
the RMF+RRPA calculations without pairing, with pairing correlations included only in the 
RHB ground state (no dynamical pairing), and with the fully self-consistent RHB+RQRPA 
response. A large configuration space enables the separation of the zero-energy mode that 
corresponds to the spurious center of mass motion. In the present calculation for 22 this 
mode is found at E = 0.04 MeV. 

The isovector dipole response in neutron-rich oxygen isotopes has recently attracted con- 
siderable interest because these nuclei might be good candidates for a possible identification 
of the low- lying collective soft mode (pygmy state), that corresponds to the oscillations of 
excess neutrons out of phase with the core composed of an equal number of protons and neu- 
trons (31L . The strength functions shown in Fig. |21 illustrate the importance of including 
pairing correlations in the calculation of the isovector dipole response. Pairing is, of course, 
particularly important for the low-lying strength below 10 MeV. The inclusion of pairing 
correlations in the full RHB+RQRPA calculation enhances the low-energy dipole strength 
near the threshold. For the main peak in the low-energy region (^ 8.65 MeV), in the right 
panel of Fig. |21 we display the proton and neutron transition densities. In contrast to the 
well known radial dependence of the IVGDR transition densities (proton and neutron densi- 
ties oscillate with opposite phases, the amplitude of the isovector transition density is much 
larger than that of the isoscalar component), the proton and neutron transition densities 
for the main low-energy peak are in phase in the nuclear interior, there is no contribution 
from the protons in the surface region, the isoscalar transition density dominates over the 
isovector one in the interior, and the strong neutron transition density displays a long tail 
in the radial coordinate. A similar behavior has been predicted for the light neutron halo 

n 

nuclei 6 He, n Li and 12 Be in Ref. |33j], where it has been shown that the long tails of the wave 
functions of the loosely-bound neutrons are responsible for the different radial dependence 
of the transition densities that correspond to the soft low-energy states as compared to those 
of the giant resonances. 

The effect of pairing correlations on the isovector dipole response in 22 is very similar to 
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the one obtained in the HFB+QRPA framework (Fig. 8 of Ref. jjj). In the low-energy region 
below 10 MeV, however, the pairing interaction used in the QRPA calculation produces a 
much stronger enhancement of the dipole strength, as compared to the results shown in 
Fig. |21 The reason probably lies in the choice of the pairing interaction. While we use the 
volume Gogny pairing, in Ref. jjj a density- dependent delta force was used in the pp channel. 
This interaction is surface peaked and therefore produces a stronger effect on the low-energy 
dipole strength near the threshold. Nevertheless, we emphasize that the RHB+RQRPA 
results for the low-lying dipole strength distribution in 22 are in very good agreement with 
recent experimental data |32| . 

In the left panel of Fig.|H]we display the RHB+RQRPA isoscalar and isovector quadrupole 
(J 71 " = 2 + ) strength distributions in 22 0. The low-lying J n = 2 + state is calculated at 
E = 2.95 MeV, and this value should be compared with the experimental excitation energy 



of the first 2+ state: 3.2 MeV 



3J]. The strong peak at E = 22.3 MeV in the isoscalar 



strength function corresponds to the isoscalar giant quadrupole resonance (IS GQR). The 
isovector response, on the other hand, is strongly fragmented over the large region of excita- 
tion energies E ~ 18 — 38 MeV. The effect of pairing correlations on the isoscalar response 
is illustrated in the right panel of Fig. El where again the full RHB+RQRPA strength func- 
tion is compared to the RMF+RRPA calculation without pairing, and with the response 
obtained when the pairing interaction is included only in the RHB ground state (no dy- 
namical pairing). As one would expect, the effect of pairing correlations is not particularly 
pronounced in the giant resonance region. The inclusion of pairing correlations, however, 
has a relatively strong effect on the low-lying 2 + state. This is seen more clearly in the left 
panel of Fig. HI where only the low-energy portion of the isoscalar strength distributions in 
22 O is shown. With respect to the RRPA calculation, the inclusion of the pairing interaction 
in the static solution for the ground state increases the excitation energy of the lowest 2 + 
state by ~ 3 MeV. The fully self-consistent RHB+RQRPA calculation lowers the excitation 
energy from w 4.5 MeV to E — 2.95 MeV. The inclusion of pairing correlations increases the 
collectivity of the low-lying 2 + state. A very similar result for the low-lying quadrupole state 
in 24 has been obtained by Matsuo in the HFB+QRPA framework P|. The proton and 
neutron transition densities for the 2 + state at E = 2.95 MeV are shown in the right panel 
of Fig. 01 They display a characteristic radial dependence. Both transition densities are, 
of course, peaked in the surface region, but the proton contribution is much smaller. The 
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RHB+RQRPA results for the 2 + excitations are in agreement with non-re 



calculations of the quadrupole response in neutron rich oxygen isotopes 



ativistic QRPA 
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IV. EVOLUTION OF THE LOW-LYING ISOVECTOR DIPOLE STRENGTH IN 
SN ISOTOPES AND N=82 ISOTONES 

The dipole response of very neutron-rich isotopes is characterized by the fragmentation 
of the strength distribution and its spreading into the low-energy region, and by the mixing 
of isoscalar and isovector modes. It appears that in most relatively light nuclei the onset of 
dipole strength in the low-energy region is due to non-resonant independent single particle 
excitations of the loosely bound neutrons. The structure of the low-lying dipole strength 
changes with mass. As we have shown in the RRPA analysis of Ref. [lOj, in heavier nuclei 
low-lying dipole states appear that are characterized by a more distributed structure of the 
RRPA amplitude. Among several peaks characterized by single particle transitions, a single 
collective dipole state is identified below 10 MeV, and its amplitude represents a coherent 
superposition of many neutron particle-hole configurations. 

Very recently experimental data have been reported on the concentration of electric dipole 
strength below the neutron separation energy in N = 82 semi-magic nuclei. The distribution 
of the electric dipole strength in 138 Ba, 140 Ce, and 144 Sm displays a resonant structure 
between 5.5 MeV and 8 MeV, exhausting « 1% of the isovector El EWSR 0. In 138 Ba 
negative parity quantum numbers have been assigned to 18 dipole excitations between 5.5 



MeV and 6.5 MeV 



38|. 



In Figs. El and El we display the isovector dipole strength distributions in eight N = 82 
isotones, calculated in the RHB+RQRPA framework with the NL3+D1S combination of 
effective interactions. The calculation is fully self-consistent, with the Gogny finite-range 
pairing included both in the RHB ground state, and in the RQRPA residual interaction. The 
isovector dipole response is shown for even-Z nuclei from 146 Gd to the doubly magic 132 Sn. 
In addition to the characteristic peak of the isovector giant dipole resonance (IVGDR) at 
« 15 MeV, the evolution of the low- lying dipole strength with decreasing proton number is 
clearly observed below 10 MeV. The strength of the low-lying dipole response increases with 
the relative increase of the neutron contribution, i.e. with reducing the number of protons. 
For the main peaks in the low-energy region below 10 MeV, in the panels on the right 
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side of Figs. 03 and |H1 we display the corresponding neutron and proton transition densities. 
The radial dependence is very different from that of the transition densities of the IVGDR 
peak. For all eight nuclei the main peak below 10 MeV does not correspond to an isovector 
excitation, i.e. the proton and neutron transition densities have the same sign. The relative 
contribution of the protons in the surface region decreases with reducing the proton number. 
In particular, for the nuclei shown in Fig.|3 138 Ba, 136 Xe, 134 Te and 132 Sn, there is practically 
no proton contribution to the transition density beyond 6 fm. The dynamics is that of a 
pygmy resonance: the neutron skin oscillates against the core. In Ref. |37( it was emphasized 
that the observed low-lying dipole states in the N=82 isotones are not just statistical El 
excitations sitting on the tail of the GDR, but represent a fundamental structure effect. In 
Fig. [7]we show that this is also the case for the RHB+RQRPA results. For the dipole strength 
distribution of 140 Ce, shown in the left panel, in the right column we compare the neutron 
and proton transition densities for the IVGDR peak at 14.31 MeV, for the peak at 12.51 
MeV, and for the main peak in the low-energy region at 8.22 MeV. The peak at 12.51 MeV, 
as well as other peaks in the interval 10-14 MeV, displays transition densities very similar to 
those of the GDR peak, i.e. these states belong to the tail of the GDR. The dynamics of the 
low-energy mode at 8.22 MeV, on the other hand, is very different: the proton and neutron 
transition densities are in phase in the nuclear interior, there is almost no contribution 
from the protons in the surface region, the isoscalar transition density dominates over the 
isovector one, and the peak of the strong neutron transition density in the surface region is 
shifted toward larger radii. 

On a quantitative level, the present RHB+RQRPA calculation does not compare too well 
with the experimental data on the low- lying dipole strength in the N=82 isotones. First, 
while the observed low-energy dipole states in 138 Ba, 140 Ce, and 144 Sm are concentrated 
between 5.5 MeV and 8 MeV, the calculated pygmy states in these nuclei are above 8 
MeV. This can be p artly explained by the low effective nucleon mass of the NL3 mean- 
field interaction [3J|. On the other hand, the excitation energies of the IV GDR are, as 
will be shown below in the example of Sn isotopes, rather well reproduced by the NL3 
interaction. The fact that NL3 reproduces the IV GDR, but not the centroid of the low- 
energy dipole strength, might indicate that the isovector channel of this force needs a better 
parameterization. Second and more important, the number of RQRPA peaks below 10 
MeV, for the operator ()59|) . is much smaller than the number of observed dipole states 
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in the low-energy region [371, 38]. The observed low-lying El strength consists of many 
states of different origin. This has been discussed in Ref. J37J]. In addition to the two- 
phonon and three-phonon states, and the soft pygmy state, in this energy region one could 
a. 8 o expect some — ional low-lying Scalar dipole strength Q, maybe mixed with 
toroidal states |41j , as well as the El strength generated by the breaking of the isospin 
symmetry due to a clustering mechanism 42]. A detailed investigation of the nature of all 
observed low-lying dipole states in N=82 nuclei is, of course, beyond the scope of the present 
analysis, since our model space does not include multi-phonon configurations. 

The Sn isotopes present another very interesting example of the evolution of the low-lying 



dipole strength with neutron number |43J. In Ref. |10| we have performed an analysis of 
the isovector dipole response of neutron-rich Sn isotopes in the relativistic RPA framework. 
The RMF+RRPA calculation has shown that, among several dipole states in the low-energy 
region between 7 MeV and 9 MeV, and characterized by single particle transitions, a single 
state is found with a more distributed structure of the RRPA amplitude, exhausting approx- 
imately 2% of the EWSR. The results of the fully self-consistent RHB+RQRPA calculation, 
with the NL3+D1S combination of effective interactions, are shown in Figs. [HI and [HI the 
isovector dipole strength functions of the Sn isotopes (left panels), and the corresponding 
proton and neutron transition densities for the main peaks in the low-energy region (right 
panels). With the increase of the number of neutrons a relatively strong peak appears below 
10 MeV, characterized by the dynamics of the pygmy resonance (see the transition densi- 
ties). The low-energy pygmy peak is most pronounced in 124 Sn. It does not become stronger 
by further increasing the neutron number, and additional fragmentation of the low-lying 
strength is observed in 132 Sn. For the Sn isotopes we can compare the RHB+RQRPA re- 
sults with experimental data on IV GDR. In the upper panel of Fig. El the experimental 
IVGDR excitation energies Q| are shown in comparison with the calculated -Egdr- The 
energy of the resonance is defined as the centroid energy 

B=^, (60) 
m 

with the energy weighted moments for discrete spectra 

m k = Y,B(J,u u )E k u . (61) 

V 

For k = 1 this equation defines the energy weighted sum rule (EWSR). The calculated 
energies of the IV GDR are in excellent agreement with experimental data, and the mass 
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dependence of the excitation energies is reproduced in detail. In the middle panel of Fig. [TO] 
we plot the calculated energies of the pygmy states. In comparison with the IV GDR, the 
excitation energies of the pygmy states decrease more steeply with the mass number. The 
ratio of the energy weighted mi moments calculated in the low (E<10 MeV) and high (E>10 
MeV) energy regions, as function of the mass number, is plotted in the lower panel of Fig. El 
The relative contribution of the low-energy region increases with the neutron excess. The 
ratio mi 5 Low/ hi^high reaches a maximum « 0.06 for 124 Sn, and it slowly decreases to ~ 0.05 
for 132 Sn. 

V. SUMMARY 

In this work we have formulated the relativistic QRPA in the canonical single-nucleon 
basis of the relativistic Hartree-Bogoliubov (RHB) model. The RHB model presents the 
relativistic extension of the Hartree-Fock-Bogoliubov framework, and it provides a unified 
description of mean-field and pairing correlations. A consistent and unified treatment of 
the ph and pp channels is very important for weakly bound nuclei far from stability. In the 
RHB framework the ground state of a nucleus can be written either in the quasiparticle 
product of independent quasi-particle states, or in the canonical basis as a highly 
correlated BCS-state. By definition, the canonical basis diagonalizes the density matrix and 
it is always localized. It describes both the bound states and the positive-energy single- 
particle continuum. The QRPA model employed in this work is fully self-consistent. For 
the interaction in the particle-hole channel effective Lagrangians with nonlinear meson self- 
interactions are used, and pairing correlations are described by the pairing part of the 
finite range Gogny interaction. Both in the ph and pp channels, the same interactions are 
used in the RHB equations that determine the canonical quasiparticle basis, and in the 
matrix equations of the RQRPA. This is very important, because the energy weighted sum 
rules are only satisfied if the pairing interaction is consistently included both in the static 
RHB and in the dynamical RQRPA calculations. The two-quasiparticle configuration space 
includes states with both nucleons in the discrete bound levels, states with one nucleon in 
the bound levels and one nucleon in the continuum, and also states with both nucleons in 
the continuum. The RQRPA configuration space includes the Dirac sea of negative energy 
states. In addition to the configurations built from two-quasiparticle states of positive 
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energy, the RQRPA configuration space contains pair-configurations formed from the fully 
or partially occupied states of positive energy and the empty negative-energy states from 
the Dirac sea. The inclusion of configurations built from occupied positive-energy states and 
empty negative-energy states is essential for the decoupling of spurious states. 

The RHB+RQRPA approach has been tested in the example of the isoscalar monopole, 
isovector dipole and isoscalar quadrupole excitations of 22 0. The NL3 parameterization has 
been used for the RMF effective Lagrangian, and the Gogny D1S finite range interaction 
has been employed in the pp channel. In the present numerical implementation the RHB 
eigenvalue equations, the Klein-Gordon equations for the meson fields, and the RQRPA ma- 
trix equations are solved by expanding the nucleon spinors and the meson fields in a basis 
of eigenfunctions of a spherical harmonic oscillator. The calculations have illustrated the 
importance of a consistent treatment of pairing correlations in the RHB+RQRPA frame- 
work. The results have been compared with calculations performed in the non-relativistic 
continuum QRPA based on the coordinate state representation of the HFB framework. It 
has been shown that the RHB+RQRPA results are in agreement with recent experimental 
data and with non-relativistic QRPA calculations of the multipole response of neutron rich 
oxygen isotopes. 

The RHB+RQRPA has been employed in the analysis of the evolution of the low-lying 
isovector dipole strength in Sn isotopes and N=82 isotones. The analysis is motivated by very 
recent data on the concentration of electric dipole strength below the neutron separation 
energy in N = 82 semi-magic nuclei. It has been shown that in neutron rich nuclei a 
relatively strong peak appears in the dipole response below 10 MeV, with a QRPA amplitude 
characterized by a coherent superposition of many neutron quasiparticle configurations. The 
dynamics of this state corresponds to that of a pygmy dipole resonance: the oscillation 
of the skin of excess neutrons against the core composed of an equal number of protons 
and neutrons. It should be emphasized that, even though the IV GDR excitation energies 
calculated with the NL3 effective interaction are in excellent agreement with experimental 
data on Sn isotopes, the pygmy peaks in the low-energy region do not compare too well with 
the data on low-lying dipole strength in N = 82 isotones. The calculated peaks are ~ 2 MeV 
higher than the experimental weighted mean energies. This might indicate that there are 
problems with the isovector channel of the effective interaction and with the effective mass. 
Namely, if the pygmy resonance is directly related to the thickness of the neutron skin, the 
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splitting between the excitation energies of the pygmy state and the IV GDR should be 
determined by the isovector channel of the effective force. A detailed quantitative analysis 
of the empirical low-lying isovector dipole response of neutron rich N = 82 nuclei in the 
RHB+RQRPA framework will be included in a forthcoming publication. 

Summarizing, the relativistic QRPA formulated in the canonical basis of the RHB model 
represents a significant contribution to the theoretical tools that can be employed in the 
description of the multipole response of unstable weakly bound nuclei far from stability. 



VI. APPENDIX: NUMERICAL DETAILS OF THE SOLUTION OF THE RQRPA 
EQUATIONS IN THE CANONICAL BASIS 

The relativistic quasi-particle RPA equations can be simplified considerably by employing 



3, 



any RHB wave 



the canonical basis. According to the theorem of Bloch and Messiah 
function can be expressed either in the quasiparticle basis as a product of independent 
quasiparticle states, or in the canonical basis as a highly correlated BCS-state. For systems 
with an even number of particles we have 

m = nK+^44)i->- (62) 

K>0 

|— ) denotes the nucleon vacuum, the operators a], and aj= create nucleons in the canonical 
basis. The occupation probabilities are given by 

^ = ;|l- £ *~ m " A ., 1- (63) 
(e„ - m - A) 2 + A 2 



e K = {K\h£,\Kj and A K = (k\A\k) are the diagonal elements of the Dirac single-particle 
Hamiltonian and the pairing field in the canonical basis, respectively. In contrast to the 
BCS framework, however, neither of these fields is diagonal in the canonical basis. The 
basis itself is specified by the requirement that it diagonalizes the single-nucleon density 
matrix p(r, r') = J2k>o Hi( r )^fc ( r ')- The transformation to the canonical basis determines 
the energies and occupation probabilities of single-nucleon states, that correspond to the 
self-consistent solution for the ground state of a nucleus. Since it diagonalizes the density 



matrix, the canonical basis is loca 
energy single-particle continuum 



ized. It describes both the bound states and the positive- 
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Many of the eigenvalues (JfiHj) of the density matrix are identically zero. In particular, 
those at very high energies in the continuum, but also those that correspond to the levels in 
the Dirac sea {no-sea approximation). Because of this degeneracy the levels in the canonical 
basis are not uniquely determined by the numerical diagonalization of the density matrix 
p(r, r'). In addition to the well defined eigenstates \k) with non-degenerate eigenvalues 
< v\ < 1, there is one set of eigenstates with eigenvalues equal and another set of 
eigenstates with eigenvalues equal 1. Any linear combination of eigenstates with eigenvalue 
(1) is again an eigenstate with eigenvalue (1). The diagonal pairing matrix elements 
A vanish in these degenerate subspaces. The corresponding single particle energies e M , 
however, are arbitrary and unphysical. Within these two subspaces the canonical basis is 
not uniquely defined. 

We therefore introduce an additional requirement, that the canonical basis in each of these 
subspaces diagonalizes the single particle Hamiltonian hp. In practical applications one thus 
first diagonalizes the matrix p. This gives all the canonical basis states with < v 2 R < 1, and 
in addition two sets of degenerate eigenstates with eigenvalues and 1. Two eigenstates \k) 
and | A) are considered degenerate if the corresponding eigenvalues differ less than a given 
parameter e^: 

\vl-v\\<e d . (64) 

In the second step the single particle Hamiltonian hp is diagonalized in the subspace of 
degenerate eigenvectors of the density matrix with eigenvalues (1). These new vectors 
are also eigenvectors of p with eigenvalues (1). This procedure uniquely determines the 
energies e K and occupation probabilities v\ of single-particle states, that correspond to the 
self-consistent solution for the ground state of a nucleus. An appropriate choice, of course, 
has to be made for the parameter e^. If it is too large, a linear combination of the eigenstates 
\k) and | A) that diagonalizes hp, will no longer be an eigenvector of the density matrix p. 

It is important to illustrate how the RQRPA results depend on the choice of the parameter 
€d in Eq. (I64|) . For the nucleus 22 0, in Fig. Illf a) we display the isovector dipole strength 
distributions, calculated with = 10~ 4 — 10~ 7 . For any two values of > 10~ 6 the 
corresponding strength distributions show pronounced differences. When < 10~ 6 , the 
dipole response does not depend any longer on its precise numerical value, and the spurious 
Nambu-Goldstone 1~ mode is found at an excitation energy < 0.1 MeV. 

The RQRPA matrix is diagonalized in the finite dimensional two-quasiparticle (2qp) 
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vector space. There are two types of 2qp states: 1) those built from qp states of positive 
energy, and 2) those formed by one fully or partially occupied state of positive energy and one 
empty negative energy state from the Dirac sea. The dimension of the RQRPA configuration 
space is thus determined by two cut-off parameters: Eq p is the maximum value of the sum 
of the diagonal matrix elements of H 11 (|54|) for the first type of 2qp states, and Ec a is the 
maximum absolute value of the sum of the diagonal matrix elements of if 11 (J54|) for 2qp 
states with one quasiparticle in the Dirac sea. The choice of the two cut-off parameters Eq p 
and E Ca is restricted by the following conditions: (a) there should be no response to the 
number operator, i.e. the Nambu-Goldstone + mode associated with the nucleon number 
conservation should have zero excitation energy, (b) the spurious excitation corresponding 
to the translation of the nucleus decouples as a zero-energy excitation mode, and (c) the 
response function does not depend on the precise numerical values of Eq p and Ec a - 

In Fig. lTlT b) we show how the response to the neutron number operator for 22 varies with 
the cut-off parameter Eq p in the range 30 - 270 MeV. The choice £^0=1700 MeV includes 
the entire negative-energy Dirac spectrum. The response is obviously reduced as the number 
of 2qp configurations increases. Already for Ec p = 90 MeV the Nambu-Goldstone + mode 
converges to < 0.1 MeV. 

A large configuration space is also necessary in order to bring the spurious 1~ state at 
zero excitation energy. In Fig. UTT c) we illustrate the convergence of the energy of the Y~ 
spurious state in 22 and 120 Sn. The excitation energies are plotted as functions of the 
energy cut-off parameter Eq p - Ec a is kept at 1700 MeV. 

The choice of the cut-off parameter Ec a has a pronounced influence on the cal- 
culated isoscalar monopole response. This is illustrated in Fig. [TlT d). where we 
show how the energies of the giant monopole resonance (GMR) in 22 and 120 Sn 
depend on the value of Ec a - For Ec a < 1150 MeV, only positive-energy 2qp 
states are included in the RQRPA basis and the excitation energies of the GMR 
peaks are simply too low. As Ec a is increased to include the negative energy 
states, the GMR excitation energies also increase and saturate for Ec a >1500 MeV. 
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FIG. 1: The strength function for the neutron number operator (left), and the isoscalar strength 
function for the monopole operator (right) in 22 O. The curves correspond to the RMF+RRPA 
calculation without pairing (dotted), with pairing correlations included in the RHB calculation of 
the ground state, but not in the RRPA residual interaction (dashed), and to the fully self-consistent 
RHB+RQRPA calculation (solid). 

FIG. 2: The isovector strength function of the dipole operator in 22 (left). The fully self-consistent 
RHB+RQRPA response (solid line) is compared with the RMF+RRPA calculation without pairing 
(dotted line), and with the RHB+RRPA calculation that includes pairing correlations only in the 
ground state (dashed line). The proton and neutron transition densities for the peak at E = 8.65 
MeV are shown in the right panel. 

FIG. 4: Low-energy portion of the isoscalar quadrupole strength distribution in 22 (left). The 
neutron and proton transition densities for the = 2 + state at E = 2.95 MeV (right). 

FIG. 5: RHB+RQRPA isovector dipole strength distributions in 146 Gd, 144 Sm, 142 Nd and 140 Ce, 
calculated with the NL3+D1S effective interaction. The corresponding proton and neutron transi- 
tion densities for the main peak in the low-energy region below 10 MeV are displayed in the panels 
on the right side. 

FIG. 6: Same as in Fig. but for the N = 82 isotones: 138 Ba, 136 Xe, 134 Te and 132 Sn. 

FIG. 7: The isovector dipole strength distribution in 140 Ce (left panel). The neutron and proton 
transition densities for the IVGDR peaks at 14.31 MeV, 12.51 MeV, and for the main peak in the 
low-energy region at 8.22 MeV (right). 

FIG. 3: The RHB+RQRPA isoscalar and isovector quadrupole strength distributions in 22 (left 
panel). In the right panel the full RHB+RQRPA isoscalar strength function (solid) is compared to 
the RMF+RRPA calculation without pairing (dotted), and with the response obtained when the 
pairing interaction is included only in the RHB ground state (dashed). 
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FIG. 8: RHB+RQRPA isovector dipole strength distributions in Sn isotopes, calculated with the 
NL3+D1S effective interaction. The corresponding proton and neutron transition densities for the 
main peak below the IVGDR are displayed in the panels on the right side. 



FIG. 9: Same as in Fig. |HJ but for the heavier Sn isotopes. 



FIG. 10: In the upper panel the experimental IV GDR excitation energies of the Sn isotopes are 
compared with the RHB+RQRPA results calculated with the NL3+D1S effective interaction. The 
calculated energies of the pygmy states are shown in the middle panel. The values of the ratio 
m i,LOw/ m i,HiGH) of the energy weighted moments mi in the low-energy region (E<10 MeV) and 
in the region of giant resonances (E>10 MeV), are plotted in the lower panel. 



FIG. 11: (a) The RQRPA isovector dipole response in 22 O calculated for different values of the 
parameter (|64|). (b) Neutron number operator response in 22 computed for four values of the 
cut-off energy parameter Eq v - (c) The position of the spurious 1~ state in 22 and 120 Sn as a 
function of the 2qp cut-off energy parameter Eq p - (d) The excitation energies of the ISGMR in 
22 O and 120 Sn as functions of the cut-off energy parameter Ec a - See text for description. 
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